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ABSTRACT 

We investigate the photometric variability of magnetized stars, particularly neu- 
tron stars, accreting through a magnetic Rayleigh-Taylor-type instability at the disk- 
magnetosphere interface, and compare it with the variability during stable accretion, 
with the goal of looking for possible quasi-periodic oscillations. The lightcurves during 
stable accretion show periodicity at the star's frequency and sometimes twice that, due 
to the presence of two funnel streams that produce antipodal hotspots near the mag- 
netic poles. On the other hand, lightcurves during unstable accretion through tongues 
penetrating the magnetosphere are more chaotic due to the stochastic behaviour of 
the tongues, and produce noisier power spectra. However, the power spectra do show 
some signs of quasi-periodic variability. Most importantly, the rotation frequency of 
the tongues and the resulting hotspots is close to the inner-disk orbital frequency, 
except in the most strongly unstable cases. There is therefore a high probability of 
observing QPOs at that frequency in longer simulations. In addition, the lightcurves 
in the unstable regime show periodicity at the star's rotation frequency in many of the 
cases investigated here, again except in the most strongly unstable cases which lack 
funnel flows and the resulting antipodal hotspots. The noisier power spectra result in 
the fractional rms amplitudes of the Fourier peaks being smaller. 

We also study in detail the effect of the misalignment angle between the rotation 
and magnetic axes of the star on the variability, and find that at misalignment angles 
> 25°, the star's period always appears in the lightcurves. 

Key words: accretion, accretion discs; instabilities; MHD; stars:neutron; stars: oscil- 
lations; stars: magnetic fields 



1 INTRODUCTION 

Accretion to magnetized stars occurs in different systems, in- 
cluding classical T Tauri stars (CTTSs), which are the progen- 
itors of solar-type stars (e.g., Bouv ier et al.|2007| , in magne- 
tized cataclysmic variables, which are accreting white dwarfs 
(e.g., |Wamer| |1995), and in millisecond pulsars, which are 
weakly magnetized accreting neutron stars (e.g., van der Klis| 
[2004) . The light curves of these stars are expected to show 
definite periodicity corresponding to the star's rotation pe- 
riod, but in addition to this they often show interesting fea- 
tures like quasi-periodic oscillations (QPOs), or lack any pe- 
riodicity altogether. The reasons for this behaviour are not 
yet perfectly clear. The lack of periodicity in many CTTSs was 
thought to be due to the lack of a dynamically important mag- 
netic field, but recent observations have shown that that is 
not always the case (e.g., Donati et al. 2007 ). QPO models in 
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neutron star systems usually attempt to look towards charac- 
teristic frequencies in the disk, and their beats with the stellar 
rotation frequency, for explanations ( van der Klis 2004). The 
variability characteristics of dwarf novae, which are a sub- 
class of cataclysmic variables, show close similarity to those 
of neutron stars (Wa rner et al.|2003] l. These observations hint 
at a common explanation for the variability features in these 
systems, perhaps based on the properties of the accretion flow 
itself. This leads us to investigate the effect of the accretion 
flow on the variability. 

Accretion to magnetized stars can be in one of two 
regimes: (i) the stable regime, in which the accretion disk is 
stopped by the star's magnetic field, and the accreting matter 
flows around the magnetosphere, forming two funnels that 
deposit matter near the star's magnetic poles (Ghosh & Lamb 
~ Konigl|[l99l] |Koldoba et al.||2002[ |Romanova et al 



1979 



2003 



2004 1, or (ii) the unstable regime, in which Rayleigh- 



Taylor-type instabilities occur at the disk-magnetosphere in- 
terface ([Arons & Lea||l976l |Elsner & Lairib|[T977l |Wang &| 
|Robertson|1985 1. The instabilities produce tongues of matter 
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that penetrate the magnetosphere and deposit matter much 
closer to the star's equator ( Ro manova et al.|2 008 ; Kulk arni &] 
|Romanova|2008| hereafter KR08) . The tongues rotate around 
the star with an angular velocity close to the inner-disk orbital 
frequency, and their number is of the order of a few. In KR08 
it was noticed that the lightcurves during unstable accretion 
can lack the clear periodicity seen in the lightcurves during 
stable accretion (see, e.g., Poutanen & Gierlihski 2003; Kulka- 
rni & Romanova 2005, hereafter KR05). Here we investigate 
this behaviour in more detail in accreting magnetized stars, 
particularly millisecond pulsars. Our goal is to look for signs 
of QPOs in the unstable regime. 

KR08 also found that the instability is suppressed at high 
misalignment angles between the star's rotation and mag- 
netic axes. Here we perform a detailed investigation of the 
boundary between stable and unstable regimes as a function 
of the misalignment angle and stellar rotation rate, focusing 
on whether the lightcurves show definite periodicity. 

We start out with a summary of our previous 3D MHD 
simulation results in Sj2] and a description of our variability 
model in ij3] This is followed in Qby a comparison of variabil- 
ity features during stable and unstable accretion, including a 
discussion of the instability on the accretion rate, misalign- 
ment angle and stellar rotation rate. We also discuss signs of 
QPO-like behaviour. We end in fj5]with our conclusions and a 
discussion of the implications of this work for observations of 
accreting millisecond pulsars. 

2 SUMMARY OF 3D MHD SIMULATION RESULTS 

For context, we summarize the model and results of our 3D 
MHD simulations described in KR08. 

2.1 Model 

The model we use is the same as in our earlier 3D MHD sim- 
ulations (Koldoba et al. 2002; Romanova et al. 2003, 2004). 
The star has a dipole magnetic field, the axis of which makes 
an angle O with the star's rotation axis. The rotation axes 
of the star and the accretion disk are aligned. There is a low- 
density corona around the star and the disk which also rotates 
about the same axis. To model stationary accretion, the disk is 
chosen to initially be in a quasi-equilibrium state, where the 
gravitational, centrifugal and pressure gradient forces are in 
balance (Romanova et al. 2002). Simulations were done for 
both relativistic and non-relativistic stars. General relativis- 
tic effects, which are important for neutron stars, are mod- 
elled using the Paczyhski-Wiita potential = GM/(r — r g ) 
(Paczyriski & Wiita 1980), where M is the mass of the star 
and r g = 2GM/c 2 is its Schwarzschild radius. This poten- 
tial reproduces some important features of the Schwarzschild 
spacetime, such as the positions of the innermost stable and 
marginally bound circular orbits. Viscosity is modelled using 
the a-model (Shakura & Sunyaev 1973; Novikov & Thorne 
1973), and controls the accretion rate through the disk. To 
model accretion, the ideal MHD equations are solved numer- 
ically in three dimensions, using a Godunov-type numerical 
code, written in a "cubed-sphere" coordinate system rotating 
with the star (Koldoba et al. 2002; Romanova et al. 2003). 
The boundary conditions at the star's surface amount to as- 
suming that the infalling matter passes through the surface of 



the star. So the dynamics of the matter after it falls on the star 
is ignored. The inward motion of the accretion disk is found to 
be stopped by the star's magnetosphere at the Alfven or mag- 
netospheric radius, where the magnetic and matter energy 
densities become equal. The subsequent evolution depends 
on the parameters of the model. 

2.2 Reference Values 

The simulations are done using dimensionless variables. For 
every physical quantity q, the dimensionless value is defined 
as q — q/qo, where qo is the reference value for q. Because of 
the use of dimensionless variables, the results are applicable 
to a wide range of objects and physical conditions, each with 
its own set of reference values, provided the magnetospheric 
radius r m is upto a few times the radius of the star R; r m is 
determined by the balance of the magnetospheric and matter 
pressures, so that the modified plasma parameter at the disk- 
magnetosphere boundary /J — (p + pv 2 )/(B 2 /8tt) ^ 1. To 
apply our simulation results to a particular situation, we have 
the freedom to choose three parameters, and all the reference 
values are calculated from those. We choose the mass, radius 
and equatorial surface magnetic field of the star as the three 
independent parameters. Appendix [X] shows how the refer- 
ence values are determined, and lists the reference values for 
three classes of central objects — classical T Tauri stars, white 
dwarfs and neutron stars. 

Subsequently, we drop the primes on the dimensionless 
variables and show dimensionless values everywhere unless 
otherwise stated. 



2.3 Results 

We found that accretion to magnetized stars with magneto- 
spheres a few times the stellar radius in size can be either 
stable or unstable. In the stable regime, the accretion disk is 
stopped by the star's magnetic field, and the accreting matter 
flows around the magnetosphere, forming two funnels that 
deposit matter near the star's magnetic poles. This produces 
two antipodal hotspots, and the lightcurves show periodic- 
ity at either the star's frequency or twice that, depending 
on the viewing geometry. In the unstable regime, Rayleigh- 
Taylor-type instabilities occur at the disk-magnetosphere in- 
terface, producing tongues of matter that penetrate the mag- 
netosphere and deposit matter much closer to the star's equa- 
tor. The tongues rotate around the star with an angular veloc- 
ity close to the inner-disk orbital frequency, and their number 
is of the order of a few. The lightcurves often lack the clear pe- 
riodicity seen during stable accretion. However, the fact that 
the tongues rotate at approximately the inner-disk orbital fre- 
quency leads to the hope of observing some periodicity or 
quasi-periodicity at that frequency. 

The instability is associated with high accretion rates, 
and coexists with funnel flows for a relatively broad range of 
accretion rates. It is expected to occur in most accreting sys- 
tems for typical values of mass, radius, surface magnetic fields 
and accretion rates. We quote here from KR08 the empiri- 
cally obtained critical accretion rate, measured at the stellar 
surface, separating the stable and unstable accretion regimes 
(see KR08 for a more detailed discussion of instability crite- 
ria): 
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CTTSs 


White dwarfs 


Neutron stars 


r m = (2-3) R 


1.7 x 10" 7 


1.1 x 10~ 7 


1.8 x 10" 8 


r m = (4-5) R 


2.1 x 1(T 8 


1.4 x icr 8 


2.2 x icr 9 



Table 1. Characteristic critical accretion rates M cr itfi yr _1 ) for 
various types of accreting systems, for two values of the magneto- 
spheric radius r m . 
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where B, R and M are the surface magnetic field, radius and 
mass of the star respectively (in dimensional units), and the 
subscript '0' denotes the reference values shown in Table |A1| 
The characteristic critical accretion rate M crit ,o depends on 
the magnetospheric radius, which is set by the dimensionless 
value of the stellar magnetic dipole moment fi and is (4-5) R 
for fx = 2 and (2-3)R for fi — 0.5. The corresponding values 
of M cr it ,o for classical T Tauri stars, white dwarfs and neutron 
stars are given in Table ^ 



3 VARIABILITY MODEL 

We calculate lightcurves from the accretion hotspots on the 
stellar surface, taking into account the general relativistic 
effects of light bending and gravitational redshift, and the 
Doppler effect, which are important for neutron stars. Light 
travel time effects are negligible and are not taken into ac- 
count ( [Poutanen & Gierli hski 2003; KR05). To calculate the 
emission from the hotspots, we integrate the emitted flux over 
the star's surface, using the method described by KR05. The 
observed flux is then given by 



F 



<IS' 



,I' E (R',\)cos\ 
7 5 (1-/9m) 5 



(2) 



Here, primed and unprimed quantities are measured in the 
star's and observer's rest frames respectively. R and r g are the 
star's radius and Schwarzschild radius respectively. For a point 
on the star's surface rotating with velocity v, f3 = v/c and 7 = 
(1 — /3 2 ) _1,/2 . Ie is the frequency-integrated specific intensity 
of the emitted radiation. Finally, A is the angle between the 
local radial direction and the direction in which a light ray 
needs to be emitted so as to reach the observer, taking light 
bending into account. 

We now assume that the entire kinetic and thermal en- 
ergy of the matter falling onto the star's surface is converted 
into blackbody radiation, which is emitted isotropically. Then 
the flux emitted from the star's surface is 



F e (R',X)=wIe(R',X), 



giving, 



F 



1 



dS' 



cos A>0 
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(4) 



Before performing the variability analysis, one more 
thing needs to be done. Fig. [T] shows the hotspot lightcurve 
from one of our simulations. Periodic variability is clearly seen 
superposed on a slowly varying background. The variation of 
this background depends on the long-term behaviour of the 
disk in our simulations. However, this is not of interest for 
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Figure 1. Hotspot lightcurve from one of our simulations (solid line), 
the slowly varying background flux (dashed line), the difference (dot- 
ted line), and the approximate steady-state flux (solid horizontal 
line) . 
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Figure 2. Stable accretion case which has the following parameters: 
viscosity parameter a = 0.02, misalignment angle = 15°, stellar 
rotation period P = 2.0, stellar magnetic moment fi = 1. The ob- 
server inclination angle is % = 15°. Top left: Constant plasma density 
surface showing matter flow around star, and magnetospheric mag- 
netic field lines. Bottom left: Northern hotspot, as seen from above 
the magnetic pole, showing the emitted flux. Right column: Hotspot 
lightcurve (normalized) and its Fourier (normalized) and wavelet 
spectra. The colours in the hotspot and wavelet plots range from 
white (low) through blue to red (high). The numbers inside the 
Fourier and wavelet plots are values of v/v,, and the percentages 
are fractional rms amplitudes. 



the short-term variability considered in the present work. We 
therefore smooth the lightcurve to obtain the background and 
subtract it from the lightcurve before analyzing the variability. 
This makes the short-term variability much clearer. In all sub- 
sequent figures, we show only the subtracted lightcurves. We 
estimate the pulsed fraction of the emission using an approx- 
imate value for the background lightcurve flux in the steady 
state. 

We show dimensionless values for all quantities in this 
work. To convert them into physical units, they need to be 
multiplied by the reference values listed in Table |A1| In par- 
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ticular, the values of time are in units of the Keplerian period 
Po at r = i? = 3/?, and not of the stellar rotation period, 
since we consider stars with different rotation rates. The ob- 
server inclination angle varies from 0° at the north rotational 
pole of the star to 90° in the equatorial plane. 

One final note: the lightcurves presented in this work 
are calculated for strongly relativistic systems like neutron 
stars, and are therefore not applicable to other types of ac- 
creting magnetized stars like classical T Tauri stars or proto- 
stars. However, the analysis of hotspot motion on the stellar 
surface in section |4.5| uses the direct results of the 3D MHD 
simulations. These simulations, as mentioned earlier, approx- 
imate general relativistic effects using the Paczyriski-Wiita po- 
tential, but simulations without taking these effects into ac- 
count show that the instability behaves in a very similar man- 
ner. We therefore expect the analysis of hotspot motion to be 
applicable to these two types of accreting systems as well. 



4 VARIABILITY IN THE STABLE AND UNSTABLE 
REGIMES 

4. 1 Stable Accretion 

In the stable regime, the funnels produce two antipodal 
hotspots near the star's magnetic poles. The funnels are 
steady, and therefore the hotspots are almost fixed on the 
star's surface ( [Romanova et al. 2004; KR05). Depending on 
the misalignment and inclination angles, either one or both 
antipodal hotspots can be seen by the observer during each 
rotation of the star. Correspondingly, the power spectrum 
shows peaks at either v#, 2v, or both. Fig. [2] shows the 
lightcurve and its Fourier and waveler] spectra for one of our 



a l_ 




(a) without relativistic effects (b) with relativistic effects 



Figure 4. Comparison of the visible portion of the neutron star with 
and without relativistic effects. 

stable accretion cases, which has the following parameters: 
misalignment angle O — 15°, viscosity parameter a = 0.02, 
stellar rotation period P — 2.0 and stellar magnetic dipole 
moment fj, — 1. The observer inclination angle is i — 15°. 
In this case, as the star rotates, only the northern hotspot is 
visible to the observer. So the hotspot lightcurve is modulated 
at the star's rotation frequency. As the wavelet shows, there is 
no drift in the pulsation frequency. The fractional rms ampli- 
tudes are small because the hotspots are near the rotational 
axis, and also due to general relativistic effects (Pout anen &| 
|Gierlinski | 2003l KR05). 

Fig. [3] shows another stable accretion case, with 6 = 30° 
and i = 45° (second column). For this accretion geometry, 
the frequencies v, and 2^» are both visible. If the inclination 
angle i — 90°, only 2v, is seen (Fig. [3] third column). 



1 Since the wavelet transform uses a window of width At centred 
at time t to calculate the time-dependent frequency spectrum of the 
lightcurve from time ti to time ti, it is necessary to have t\ + At/2 < 



t < t2 — At/2. This results in a portion of the wavelet plot being 
cut off. The width At is inversely related to the frequency, so lower 
frequencies are more strongly affected by this restriction. 
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(a) northern hotspot 




(b) southern hotspot 



Figure 5. Hotspots around the north and south magnetic poles for 
the unstable case in Fig. [6] (discussed in j ]4.2[ l . 
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Figure 6. Similar plot to Fig. [2] for an unstable case, which has the 
following parameters: a = 0.02, = 5°, P = 2.0, n = 1. The 
observer inclination angle is i = 45°. The matter flow in the top 
left panel is shown by density contours on two vertical slices passing 
through the origin. 
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Figure 7. Similar plot to Fig.[2]for a strongly unstable accretion case, 
which has the following parameters: a = 0.2, © = 5°, P = 2.8, 
{i = 2. The observer inclination angle is i = 45°. 



hotspots still produces a peak at the star's rotation frequency 
in the lightcurve, but it is not as steady as in the stable case. 
The power spectrum shows additional peaks produced by the 
tongue component of the hotspots. These frequencies appear 
only sporadically, though, as the wavelet shows. This is be- 
cause of the stochastic nature of the tongues. 

Deeper in the unstable regime (Fig. [7}, the hotspots are 
created entirely by tongues. Therefore, their shape, intensity, 
number and position changes on the inner-disk dynamical 
timescale. The star's rotation frequency is absent from the 
power spectrum, and the Fourier and wavelet spectra are very 
chaotic. An important exception to this are cases with small 
magnetospheres, which show clear QPO features (Romanova 
|&Kulkarni|2008l ). 



One thing that needs to be noted, though, is that even 
in this relatively simple accretion geometry, it is not always 
straightforward to guess if the power spectrum will have 
peaks at or 2v, . Firstly, the hot spots are not exactly at the 
magnetic poles; rather, they are banana-shaped and surround 
the pole. Secondly, in neutron stars, due to gravitational bend- 
ing of light, more than half the star's surface is visible, caus- 
ing the hotspot on the far side of the star to be visible either 
for some or all stellar rotation phases (Fig. [4j see also |Be-| 
|loborodov|2002||Poutanen & Gierliriski|2 003 ; KR05). Thirdly, 
the two antipodal hotspots are not always identical (Fig. [5j 
identical spots would have appeared as mirror images of each 
other here), in which case the fundamental frequency is v, 
and not 2u t , even at high observer inclination angles. 

4.2 Unstable Accretion 

As we move into the unstable regime, we start seeing ac- 
cretion through both funnels and tongues. The hotspots be- 
come more chaotic (Fig. |6j. The funnel component of the 



4.3 Dependence of the Variability on the Accretion Rate 

The accretion rate is one of the most important factors deter- 
mining the accretion stability — high accretion rates favour 
the instability l |Kulkarni & Roman ova 2008). This can be 
crudely understood as follows (for a more detailed discus- 
sion, see KR08 and, e.g., Spruit et al.|1995| and |Li & Narayan| 
|2004| for instability criteria) . The matter at the inner disk 
boundary is strongly braked by the stellar magnetic field and 
is brought closer to corotation with the star. The effective 
gravitational acceleration, which is the difference between 
the gravitational and centrifugal accelerations, is one of the 
important parameters which decide if the accretion is stable 
— a stronger effective gravitational acceleration is favourable 
for the Rayleigh-Taylor instability. This acceleration is deter- 
mined by how far the inner-disk angular velocity is from be- 
ing Keplerian. When the accretion rate increases, the inner 
edge of the accretion disk moves inward, and the inner-disk 
matter becomes more strongly sub-Keplerian, increasing the 
effective gravitational acceleration and causing the instability 
to become stronger. 
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Figure 8. Equatorial (top row) and vertical (bottom row) slices showing the density distribution for various values of a. The colours range from 
white (low) through blue to red (high) . The black line denotes where the matter and magnetic energy densities are roughly equal. 



The most obvious way of controlling the accretion rate 
is through the disk density, but in our simulations, it is more 
convenient to vary the viscosity parameter a instead. We keep 
the dimensionless density fixed at pa = 1 in all our sim- 
ulations. Increasing a increases the accretion rate, bringing 
the inner edge of the accretion disk inwards, as reflected 
by compression of the magnetosphere (Fig. [8] bottom row). 
However, increasing the density also makes the instability 
stronger (KR08, §3.1). Thus, the reason the instability be- 
comes stronger with increasing a is not the higher viscosity 
itself, but the higher accretion rate and the smaller resultant 
magnetospheric radius. 

Therefore, to investigate the dependence of the instabil- 
ity and variability on the accretion rate, we vary a, keeping 
the other parameters fixed at the following values: misalign- 
ment angle O = 5°, stellar rotation period P = 2.8 and stellar 
magnetic moment fi — 2. Fig.[9]shows the wavelet and Fourier 
spectra for these cases. At a = 0.02 the accretion is stable, and 
the lightcurve has a peak at the star's rotation frequency (the 
peak at half that frequency is due to wandering of the hotspot 
itself, and is a separate issue). As we move on to a — 0.03, we 
start getting closer to the unstable regime. The wavelet spec- 
trum shows the peak at the star's rotation frequency begin- 
ning to waver. This wavering is due to slight wandering of the 
funnel stream induced by the instability, although the insta- 
bility is not strong enough to produce tongues at this stage. At 
a = 0.04 the accretion is unstable, and the peak near the star's 
rotation frequency is not steady; it is seen to drift between 
v/v, = 0.7 and v /v* — 1.3. The Fourier spectrum still shows 
a peak at v/v, = 1.0, but it is much broader and the spectrum 
has much more noise than in the stable accretion cases. No- 
tice that the hotspot in this case has a banana-shaped portion 
produced by the funnels, and features protruding from the ba- 
nana, which are created by tongues. The tongues, and hence 
these protrusions, rotate faster than the star, producing the 
high-frequency peaks at v/v, — 2.5 — 2.7 seen in the wavelet. 



There are two prominent protrusions, and they survive for a 
long time, rotating with an almost constant frequency. This is 
therefore a strong candidate for obtaining QPOs. We discuss 
this in more detail in §4.5| 

As we move to higher accretion rates, we see that the 
noise in the Fourier spectra continues to increase, and hints 
of peaks at other frequencies start appearing. Note, however, 
that even though these cases are unstable, there is a reason- 
ably clear peak at the star's rotation frequency, since accre- 
tion occurs through both funnels and tongues. It is only when 
we reach the strongly unstable case with a — 0.2 that the 
accretion occurs solely through tongues, and the star's rota- 
tion frequency disappears from the Fourier spectrum, which 
is dominated by other frequencies. At a = 0.3, it is difficult 
to discern any peaks in the Fourier spectrum. Observationally 
this would correspond to lack of pulsations. We discuss the 
implications of this in !j5] 



4.4 More Sample Cases 

Fig. [To] shows some more cases with interesting spectral be- 
haviour, for misalignment angle = 5°. The important thing 
to keep in mind here is that one of the important factors 
determing whether the instability exists is the difference be- 
tween stellar rotation frequency and the Keplerian frequency 
at the inner disk radius, as mentioned in the previous section. 
Thus, increasing the rotation period is favourable for the in- 
stability. 

The most interesting of the cases in Fig.[l0]is the one in 
the top row, which has the same parameters as in our main 
case (Fig. [9] a = 0.1), except for a very small magnetospheric 
radius (r m /R = 2). The accretion is strongly unstable, occur- 
ring solely through tongues, and since the inner-disk orbital 
frequency is higher, the tongues rotate at a much higher fre- 
quency than in the main case. This produces peaks at very 
high frequencies in the power spectrum. The interesting thing 
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Figure 9. Wavelet and Fourier spectra and hotspots for cases that differ only in the accretion rate, controlled by the a-viscosity parameter. In all 
these cases, = 5°, P = 2.8 and fi = 2. The observer inclination angle is i = 45°. The numbers in the wavelet and Fourier plots are values of 
v/u*. The colours in the hotspot and wavelet plots range from white (low) through blue to red (high). 
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Figure 10. Wavelet and Fourier spectra and hotspots for various cases, showing QPO-like features. The observer inclination angle is i = 45°. The 
numbers in the wavelet and Fourier plots are values of . The hotspots are all shown pole-on, except for the one in the top row, where the line 
of sight is in the plane of the rotational equator. The colours in the hotspot and wavelet plots range from white (low) through blue to red (high) . 



about unstable cases with weak magnetic fields is that the 
behaviour of the tongues is noticeably different. There are 
usually two tongues, and they rotate almost rigidly with an 
almost constant angular velocity (Fig. |lip , producing clear 
and steady QPOs. It must be pointed out, however, that the 



matter in the tongues does not rotate rigidly — it is only the 
pattern of gas flow that does. We discuss accretion to small 
magnetospheres elsewhere (Romanova & Kulkarni 2008). 

The second row of Fig. [To] shows another unstable case 
with two tongues that persist for a long time, rotating faster 
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Figure 11. Motion of the tongues for a case with fj, = 0.2, which has 
a small magnetosphere. Constant density surfaces and magnetic field 
lines are shown. 

than the star with an almost constant angular velocity, like 
in the case with a = 0.04 described in §4.3| This is another 
strong candidate for obtaining QPOs. 

The third row of Fig. [To] shows an interesting phe- 
nomenon — the instability exists only between t = 7 and 
t — 14. It is only during this period that the star's rotation fre- 
quency is absent from the wavelet spectrum. This may be rel- 
evant for the intermittency of pulsations from some accreting 
millisecond pulsars: changes in the accretion rate can cause 
the accretion to switch between being stable and unstable, 
causing the star's rotation frequency to alternately appear and 
disappear from the lightcurve power spectra. We discuss this 
in© 

The wavelet plot in the fifth row shows an interesting 
feature: the existence of two relatively long-lived frequencies 
in addition to the stellar frequency. If they remain steady over 
longer periods of time, they can be expected to produce quasi- 
periodic oscillations. 

The sixth and seventh rows show that even for relatively 
low values of the accretion rate, determined by a, if the star 
rotates slowly, the instability can be strong enough that the 
star's rotation period is completely absent from the lightcurve 
power spectrum. 

4.5 Search for Explanation of Frequencies and Possible 
QPOs 

We try to find an explanation for the various frequencies that 
appear during unstable accretion shown in the preceding sec- 
tions, by tracking the motion of the hotspots on the star's sur- 
face. The hotspots produced by the tongues rotate approxi- 
mately around the star's magnetic pole, so we choose a ring 
on the star's surface centered at the magnetic pole which con- 
tains the hotspots (Fig. |12p . Then we integrate over the lati- 
tude 8 the flux emitted from this ring, at each azimuthal posi- 
tion 4> along the ring and each instant of time t, in a reference 
frame rotating with the star. This integrated flux is shown in 
Fig. [13] for the cases with different values of the a-viscosity 




Figure 12. Setup for tracking the hotspot motion. The flux emitted 
from within the ring (bounded by the two circles and centered at 
the magnetic pole /^) is integrated over the latitude 9, at each lon- 
gitude and each instant of time, giving the azimuthal position of 
each hotspot as a function of time. Zero longitude is defined as the 
direction away from the rotation axis O. 



shown in Fig. [9] The streaks marked by solid lines in this plot 
are produced by the hotspots, and the slopes of the streaks 
give the angular frequency of rotation of the spots. We there- 
fore refer to these plots as "spot-omega" plots. The dashed 
lines indicate the continuation of the solid lines from = 
when they reach <j> = 360°. The pair of numbers adjacent 
to each solid line denote the spot rotation frequency in units 
of the stellar rotation frequency (i.e., v/v*), and the num- 
ber of rotations performed by the spot during its lifetime, as 
seen by an external observer. The cases in which funnels are 
present show horizontal bands that simply indicate that the 
hotspots produced by the funnels stay at a fixed location on 
the star's surface, due to which they have v/v* = 1 as seen by 
an external observer. These spots generally last the entire du- 
ration of the simulations; the number of rotations performed 
by them is therefore not shown. The other spots are produced 
by tongues, and rotate at approximately the orbital frequency 
of the inner-disk matter, which in turn depends on the Ke- 
plerian frequency at the inner disk radius. The angular fre- 
quencies of the spots may therefore be expected to change on 
long timescales if there are secular changes in the inner disk 
radius. 

There is some agreement between the frequencies ob- 
tained from the spot-omega plots and those from the fourier 
and wavelet plots. For example, the spot-omega plot for a = 
0.04 shows vjv, w 1.3, which is close to some of the fre- 
quencies in the corresponding wavelet plot, and is also seen 
in the fourier plot, albeit not very distinctly. The spot-omega 
plot for the a — 0.08 case shows a range of frequencies be- 
tween 1.2 and 1.6, which is consistent with the peaks seen in 
the fourier plot at that range of frequencies. The agreement, 
however, is not very encouraging. There are several possible 
reasons for this: (1) The hotspots produced by the tongues 
appear and disappear sporadically, and are relatively short- 
lived; as the spot-omega plots show, each spot only completes 
of the order of a few rotations around the star (as seen by an 
external observer). (2) Different spots separated either spa- 
tially or temporally are not perfectly correlated with regard 
to their location, again as shown by the spot-omega plots; 
the individual streaks in those plots show little correlation 
in location and duration. (3) The shape and brightness of 
the spots constantly changes. This is also seen from the spot- 
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Figure 13. "Spot-omega" plots to track hotspot motion. The normalized hotspot flux as a function of the magnetic longitude <j> and time is shown 
for the cases with different a in Fig.[9] which have = 5°, P = 2.8 and = 2. The first of each pair of numbers inside the plots is the value 
of v /u, for the hotspot, and the second is the number of rotations performed by the spot during its lifetime, as seen by an external observer. The 
colours range from white (low) through blue to red (high) . 



omega plots; the width and brightness of each streak varies 
as we move along the streak. (4) The angular velocity of the 
tongues, and therefore of the spots, is not constant; differ- 
ent tongues in the spot-omega plots have different angular 
velocities, and occasionally we even see a sudden jump in 
the frequency (e.g., for a = 0.3 at t = 12.5). Due to these 
factors, the spots do not produce a coherent modulation of 
the observed flux. In fact, lack of coherent modulation, which 
consists of effects like rapid amplitude and phase changes, or 
"jumps" in the lightcurve, could produce short-duration peaks 
in the wavelet, which would appear in the Fourier spectra if 



the duration of the lightcurves is short, as in our simulations 
(~ 10 stellar rotation periods). 

However, we believe that there is still a lot of hope of 
seeing the spot-omega frequencies in the lightcurves. The 
main culprit for the absence of these frequencies from the 
fourier and wavelet plots seems to be the short duration of 
our lightcurves combined with the lack of coherent modula- 
tion mentioned above. Since the tongues have a finite life- 
time, their contribution to the lightcurve consists of wave- 
trains with random phase, amplitude and duration, and it 
is difficult to extract the frequencies of these wavetrains if 
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the lightcurves are a few tens of stellar rotation periods long. 
However, the amplitude and duration distribution of these 
wavetrains is not very broad. Experiments with much longer, 
artificially generated lightcurves consisting of such wavetrains 
show clear presence of QPOs (Bachetti & et al. 2009), due 
to the fact that in the fourier domain, the peaks due to def- 
inite signals grow faster than the noise when the length of 
the lightcurve increases. One encouraging fact about the spot- 
omega plots is that although the individual tongues are short- 
lived, the tongue frequency is relatively well-defined for a 
given set of parameters; the spread is generally not very large, 
except in the most strongly unstable cases (a > 0.2). This 
strongly increases the chances of obtaining QPOs from longer 
lightcurves. The spread in the hotspot frequency may con- 
tribute to the width of the QPO peak. The center of the peak 
would be determined by the center of this frequency distri- 
bution, which, as mentioned above, is very close to the inner 
disk frequency. We may therefore expect QPOs close to the 
inner disk frequency, except perhaps in the most strongly un- 
stable cases. 

Another interesting question concerns the significance of 
the frequencies in the fourier and wavelet plots. As mentioned 
earlier, a significant portion of them probably amount to noth- 
ing more than noise. However, it is possible that the hotspot 
shape and brightness changes have some quasi-periodicity 
that is difficult to detect using the spot-omega plots. The spot- 
omega plots do show some signs of this: the tongue spots are 
usually brighter around 4> = 0, that is, when the hotspots, ro- 
tating around the magnetic pole, come close to the disk. This 
might produce some of the peaks in the wavelet and fourier 
plots. The other peaks would be expected to contribute to 
broadband noise in the power spectrum, the effect of which 
would a smaller relative amplitude of the peak at the star's 
rotation frequency than in the stable case. It is difficult to as- 
certain this without longer lightcurves. 

It is also very interesting to observe the effect of increas- 
ing accretion rate (via increasing a) on the tongue frequen- 
cies in the spot-omega plots. The most important effect is that 
the frequencies show an increasing trend. This is due to the 
fact that the disk comes in closer and the inner-disk orbital 
frequency, which determines the tongue frequency, increases. 
We also see signs of the tongue behaviour becoming more 
chaotic; the tongue lifetimes become shorter (although this is 
not a very clear trend), and the frequency spread increases. 
These two effects would be responsible for decreasing the 
strength of the QPO peak. Thus, is reasonable to expect that 
the most strongly unstable cases (a > 0.2) would not show 
any QPOs. The lightcurve power spectra in these cases would 
completely lack pulsations. 

For stars with much smaller magnetospheres than those 
considered here, we find a different regime of accretion which 
we call the "magnetic boundary layer regime," in which the 
behaviour of the tongues is much more coherent, producing 
well-defined peaks in the Fourier and wavelet spectra. We dis- 
cuss this regime in detail in Romanova & Kulk arni| ( |2008| ). 

4.6 Dependence of the Instability on the Misalignment 
Angle 

So far we have considered cases with a small misalignment 
angle, O = 5°. For higher O the accretion geometry is slightly 
different even in the stable accretion regime, due to the 




Figure 14. Left panel: Cutaway view of the accretion flow in an un- 
stable case with = 30° . Right panel: Accretion flow in an unstable 
case with = 60°. 



stronger non-axisymmetry. The funnel flows are closer to the 
equatorial plane and produce elongated, rather than banana- 
shaped, hotspots (Roma nova et al.||2004[ KR05). It is there- 
fore interesting to explore the effect of O on the unstable ac- 
cretion geometry. At high misalignment angles (O > 25°), 
the funnels are always present. In this situation, in the un- 
stable cases, the tongues have low density, and are restricted 
to being under the edges of the funnels, due to which the 
funnels appear to be "connected", as Fig.[l4]shows. This con- 
nection appears and disappears stochastically. As we attempt 
to move into the strongly unstable regime (by either increas- 
ing the accretion rate or decreasing the star's rotation rate), 
the connecting tongues increase in density and become per- 
manent. The magnetic field lines near the edges of the funnels 
are pried apart by these tongues and stay in that configura- 
tion, creating a new equilibrium configuration with a stable 
accretion flow, albeit partially inside the magnetosphere. 

To explore the effect of on the variability properties, 
we first vary keeping the other parameters fixed at the fol- 
lowing values: viscosity parameter a — 0.1, stellar rotation 
period P = 2.8 and stellar magnetic moment fi — 2. Fig. 
|15| shows the effect on the power spectrum of increasing 6, 
which is similar to decreasing the accretion rate. For the pa- 
rameters chosen here, funnels are present for all values of 6. 
For small 6, the peak at the star's rotation frequency is visi- 
ble, but the power spectrum is noisy due to the presence of 
the tongues. As we go to higher 0, the motion of the instabil- 
ity tongues is azimuthally restricted as mentioned above. The 
hotspots are therefore more or less fixed on the star's surface, 
and so the stellar rotation frequency becomes clearer in the 
lightcurve. 

Notice the high harmonic content of the power spectrum 
for = 60°. This is probably due to the fact since the mis- 
alignment angle is large, the hotspot (which is close to the 
magnetic pole) is close to the rotational equator, and there- 
fore has a higher linear velocity, as seen by an external ob- 
server, than for smaller misalignment angles. The distortion 
of the lightcurve away from a sinusoidal shape due to the 
Doppler effect is therefore stronger (see, e.g., Pout anen &| 
|Gierliriski|2003l KR05). 

From the point of view of observations, it is interesting 
to know whether or not the lightcurves show periodicity at 
the star's rotation frequency. In order to examine in detail the 
parameter ranges for which they do, we performed a range of 
simulation runs for different values of the misalignment an- 
gle and the stellar rotation period P. We chose the stellar 
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Figure 15. Wavelet and Fourier spectra and hotspots for cases that differ only in the misalignment angle ©. In all these cases, a = 0.1, P = 2.8 
and /.t = 2. The observer inclination angle is i = 45°. The numbers in the wavelet and Fourier plots are values of u/u t . The colours in the hotspot 
and wavelet plots range from white (low) through blue to red (high) . 



magnetic moment to be fi — 1 and the viscosity parameter 
to be a — 0.02 in all these runs. Fig. |16| shows the results. 
The lightcurve lacks periodicity at the star's rotation period 
only below a critical misalignment angle Q cr it which depends 
on the rotation rate. We see that for P < 2.2, @ cr it < 0; 
the star's period is visible for all misalignment angles. As we 
decrease the star's rotation rate, the difference between the 
star's rotation rate and the Keplerian frequency at the inner 
disk radius increases, causing stronger magnetic braking of 
the matter. Hence, the effective gravitational acceleration at 
the inner disk edge increases, causing the instability to be- 
come stronger, and Q c nt to increase. However, when O c ,it 
reaches ~ 25°, it stops rising. For larger misalignment angles 
than this, the funnels are always present, and their contribu- 
tion to the hotspots gives rise to variability at the star's rota- 
tion frequency. 



4.7 Dependence of Results on Model Parameters 

One question that remains concerns the dependence of the 
results described above on the viscosity parameter a and the 
grid resolution. The effect of these parameters on the be- 
haviour of the tongues has been discused in detail in KR08. 
To summarize, one of the important factors that determine 
the appearance of the instability is the balance between the 
gravitational and centrifugal forces at the inner disk bound- 
ary, which in turn depends on the accretion rate and the stel- 
lar rotation rate and magnetic field. The viscosity parameter 
comes into play only insofar as it determines the accretion 
rate, and does not impact the variability properties directly. 
To check this, we chose a case with a small value of a where 
the accretion was stable, and increased the density in the disk, 
to find that the accretion became unstable (KR08, §3.1). Sim- 
ilarly, from the last two rows of Fig. [To] we see that in spite 
of the relatively low values of a, the instability exists and is 
so strong that there is no accretion through funnels, as evi- 
denced by the absence of the star's rotation frequency from 
the power spectra. Thus, it is the accretion rate, and not the 
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Figure 16. The regimes in which the lightcurve shows the star's rota- 
tion period, as a function of the misalignment angle (in degrees) 
and rotation period P, for a = 0.02 and fi = 1. The circles and trian- 
gles represent simulation runs that have lightcurves with and without 
the star's rotation period respectively. 



viscosity, which is important for determining the properties of 
the instability. The viscosity coefficient does, however, impact 
the radial velocity of the accreting matter; the effect, if any, of 
this on the instability needs to be explored. As far as the de- 
pendence on the grid resolution is concerned, the important 
thing to note is that the aspects of the tongue dynamics which 
determine the variability properties are the number and rota- 
tion rate of the tongues, and to a lesser degree, the location 
on the star where the tongues deposit matter. These proper- 
ties are found to be independent of the grid resolution (KR08, 
§3.4). 



5 CONCLUSIONS AND DISCUSSION 

The main results of this work are the following: 

(i) For small misalignment angles between the rotation 
and magnetic axes of the star, the tongues that penetrate the 
magnetosphere during unstable accretion, and the resulting 
hotspots on the star's surface, have an orbital frequency close 
to that of the inner-disk in most cases. We therefore expect to 
obtain quasi-periodic oscillations (QPOs) at the inner-disk or- 
bital frequency in a wide range of unstably accreting systems. 
As we move deeper into the unstable regime, the tongue be- 
haviour becomes more and more chaotic, and we may expect 
the QPOs to become weaker and gradually disappear. 

(ii) There is significant accretion through antipodal funnels 
even during unstable accretion for a broad range of parame- 
ter values. The lightcurves therefore show periodicity at the 
stellar rotation frequency or twice that, depending on the mis- 
alignment angle and the viewing geometry, as in cases with 
stable accretion. 

(iii) For large misalignment angles, the accretion funnels 
are always present, and the tongue motion is restricted, due 
to which the power spectra always show the star's rotation 
frequency. 

• At small misalignment angles (© < 25°) and mag- 
netospheric sizes of a few stellar radii, the hotspots due to 
the tongues exhibit stochastic behaviour, and thus contribute 
only noise to the power spectra on the timescales explored 
here. The important thing, however, is that the rotation fre- 
quency of the tongues and the resulting hotspots stays rel- 



atively steady at a value close to the inner-disk orbital fre- 
quency, even in the nonlinear regime of the instability, ex- 
cept in the most strongly unstable cases. This could produce 
QPOs in the power spectra at frequencies close to the inner- 
disk orbital frequency. Longer lightcurves are needed to verify 
this. The stochastic behaviour of the tongues may contribute 
to low-frequency broadband noise in the power spectra (see, 
e.g., |Burderi et al.|[1997| . We see some signs of this in our 
power spectra; the noise amplitude is lower at high frequen- 
cies in most cases. At the same time, due to the presence of 
the funnels, the star's frequency is expected to appear in the 
power spectra in most unstable cases. 

There exists a strongly unstable regime of accretion, how- 
ever, in which the funnels are absent, and the lightcurves do 
not show periodicity at the star's rotation frequency; more- 
over, the tongue behaviour is more chaotic, making the pres- 
ence of QPOs much less likely. These cases may correspond to 
accreting neutron stars without any pulsations. Since unsta- 
ble accretion occurs at relatively high accretion rates (j ]4.3| 
also KR08), the above results have a few implications: (1) 
If the accretion rate is close to the boundary between sta- 
ble and unstable regimes, slight changes in the accretion rate 
can cause the accretion to episodically switch between being 
stable and unstable, causing corresponding appearance and 
disappearance of pulsations (j ]4.4p . This might be a possible 
explanation for the behaviour of intermittent pulsars (e.g., 
|Altamirano et al.||2008l |Casella et al.| [2008). One of the at- 
tractive features of this idea comes from the fact that dur- 
ing stable accretion, the azimuthal location of funnels with 
respect to the star is fixed, even across periods of unstable 
accretion. The pulsations in intermittent pulsars would there- 
fore be expected to be coherent in phase across periods of lack 
of pulsations. This "phase memory" has been observed in the 
intermittent pulsars mentioned above. (2) Lack of pulsations 
does not imply a dynamically important magnetic field. The 
most strongly unstable cases shown in this work lack pulsa- 
tions in their lightcurves, but have distinct magnetospheres 
that strongly affect the matter flow around the star. The mag- 
netic moment in these cases is stronger by a factor of ~ 20 
than in cases that show field burial, discussed in Romanova & 
Kulkarni ( 2008 ). Field burial is an argument that is often ad- 
vanced as an explanation for the paucity of accreting millisec- 
ond pulsars. The unstable regime complements this argument 
by providing another way for pulsations to be absent. 

• At large misalignment angles (0 > 25°), funnels are 
always present, and the motion of the tongues is azimuthally 
restricted to being close to the edges of the funnels. The 
hotspots are therefore approximately fixed on the star's sur- 
face, due to which the lightcurves always show the stellar ro- 
tation frequency, even during strongly unstable accretion. 

Although this work focuses on neutron stars, we expect 
the above conclusions to be applicable to other types of ac- 
creting magnetized stars in a broad sense. One thing needs 
to be noted, however: the magnetospheres in all the simula- 
tions presented here are a few stellar radii in size. For much 
smaller magnetospheres, the variability is very different l[Ro-| 



|manova & Kulkarni 2008). In the stable regime, the funnels 
show significant wandering, due to which the star's rotation 
frequency is absent from the lightcurve power spectra. In the 
unstable regime, the tongues are much more coherent, pro- 
ducing distinct peaks in the power spectra. 
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Table Al. Sample reference values of the dynamical quantities used 
in our simulations. 



reference velocity is the Keplerian velocity at Ro, vo = 
(GM/R ) 1/2 , and Wo = vo/Ro is the reference angular ve- 
locity. The reference time is the Keplerian rotation period 
at Ro, Po = 2-kRq/vo- The reference surface magnetic field 
of the star at the magnetic equator is B* . The reference 
magnetic field, Bo, is the initial magnetic field strength at 
r = Ro, assuming a surface magnetic field of B* . The ref- 
erence magnetic dipole moment is /i = B t0 Rt = B Ro- 
A dimensionless magnetic moment of // then corresponds 
to a surface magnetic field of B* = n'B,o- The reference 
density is taken to be po = Bo/vq. The reference pressure 
is po — PoVq. The reference temperature is To = po/TZpo, 
where TZ is the gas constant. The reference accretion rate is 
Mo = PovoRo- The reference energy flux is Eo = pov'oRo- The 
reference value for the effective blackbody temperature of the 
hot spots is (T c ff)o = (po^o/"") 1 ^ 4 , where a is the Stefan- 
Boltzmann constant. Table |Al| shows sample reference values 
for three classes of objects: classical T Tauri stars (CTTSs), 
white dwarfs and neutron stars. 



APPENDIX A: REFERENCE VALUES 

As stated in section j ]2.2| our simulations are done using di- 
mensionless variables, obtained by dividing the dimensional 
variables by their respective reference values. The reference 
values are determined as follows: The unit of distance Ro 
is chosen such that the star has radius R — 0.35i?o- The 
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